Expectation does not exist: Cauchy Distribution

Make a sample from a Cauchy distribution and compute the sample mean and median.

z<- rcauchy(1000)
mean(z)

## [1] 0.5208591

summary(z)

##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## -198.2917   -0.8392   -0.0202    0.5209    0.9410  324.9145

y<- rnorm(1000,3,.5)
mean(y)

## [1] 3.010735

summary(y)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   1.250   2.652   3.019   3.011   3.348   4.420

Make many means and medians.

x.1 <- replicate(1000,mean(rcauchy(10000)))
x.2 <- replicate(1000, median(rcauchy(10000)))

summary(x.1) #summary of mean

##       Min.    1st Qu.     Median       Mean    3rd Qu.       Max. 
## -1925.7108    -0.8812     0.0439    -0.3628     1.0505  1953.1474

summary(x.2) #summary of median

##       Min.    1st Qu.     Median       Mean    3rd Qu.       Max. 
## -0.0522713 -0.0107851 -0.0001419  0.0003810  0.0117740  0.0623587

Maybe the sample size is too small for computing the mean ?

Plot the ecdf of a suitable sample.

plot.ecdf(rcauchy(1000))

Standard deviation of a Cauchy rv DOES NOT EXIST

boxplot(replicate(1000, sd(rcauchy(100))))

Example1: Second moment of a binomial random variable

n = 30
p = .4
sum((0:n)^2*dbinom(0:n,n,p)) # E(X^2) = sum(x^2 * P(X=x))

## [1] 151.2

mean(rbinom(10000,n,p)^2)  # E(X^2)

## [1] 151.0061

Example 2: Third central moment of a Poisson random variable

lambda = 4
n = 50
my.mean <- sum((0:n)*dpois(0:n, lambda)) # mu =E(X)=sum(x*P(X=x))

my.mean

## [1] 4

sum(((0:n) - my.mean)^3*dpois(0:n, lambda)) #E((X-mu)^3) = sum((x-mu)^3 * P(X=x))

## [1] 4

mean((rpois(100000,lambda) - my.mean)^3) #E((X-mu)^3)

## [1] 3.80861

Example 3: Second central moment

1. Standard deviation of a binomial distribution

n = 30
p = .4
sqrt(n*p*(1-p))

## [1] 2.683282

sd(rbinom(10000,n,p))

## [1] 2.687843

2. Standard deviation of a normal distribution

my.mu = 1
my.sigma = 3.5
sd(rnorm(1000, my.mu, my.sigma))

## [1] 3.465753

sqrt(mean((rnorm(1000,my.mu, my.sigma) - my.mu)^2)) #E((X-mu)^2))

## [1] 3.522004

Example 4: Moments package

library(moments)
x <- rnorm(100)

## Compute the mean
moment(x) # E(X)

## [1] 0.1765234

## Compute the 2nd centeral moment (= var)
moment(x, order=2,central=TRUE)

## [1] 0.8559815

## Compute the 3rd absolute centeral moment
moment(x, order=3, central=TRUE, absolute=TRUE)

## [1] 1.260265

Example 5: Mixed Distributions

Let \(X.0 \sim U(0,3)\), \(X.1 \sim Poisson(4)\), and let \(X = X.0 + X.1\)

1. Make an empirical cdf of \(X\)

2. Make an empirical cdf of \(X|X < 8\) in the same plot

3. Make an empirical cdf of \(X.1|X<8\)

N = 10000
X.0 <- runif(N,0,3)
X.1 <- rpois(N,4)
X <- X.0 + X.1

plot.ecdf(X)
X.2 <- X[X<8]
mycdf = ecdf(X.2)
t <- seq(0,10,by = .01)
lines(t,mycdf(t), col = 2)

plot.ecdf(X.1[X<8])

Example 6: Joint Distribution

\(X1 > 0\), \(X2 > 0\), \(X1 + X2 < 1\), uniform density = 2

N = 10000  # number of initial simulation
set.seed(101)
mydf.0 <- data.frame(X1 = runif(N), X2 = runif(N))
mydf.0 <- mydf.0[rowSums(mydf.0) < 1,]
head(mydf.0)

##            X1         X2
## 2  0.04382482 0.55931317
## 5  0.24985572 0.37856509
## 8  0.33346714 0.52437539
## 12 0.70687474 0.02394005
## 15 0.45512059 0.06166735
## 19 0.41166683 0.15480003

plot(mydf.0, pch = 46)

Example 7: to illustrate Law of Total Probability

1. First simulate large number of conditions \(X \in A_j\). These occur as in a geometric distribution with p = 0.5.

_Because geometric p.m.f \(=p(1-p)^{x-1}=(1/2)(1/2)^{j-1}=(1/2)^j\)

n = 100000
A.1 <- rgeom(n,.5)
table(A.1)

## A.1
##     0     1     2     3     4     5     6     7     8     9    10    11    12 
## 49857 25097 12570  6243  3165  1507   800   378   202    78    48    28    12 
##    13    14    15    16    17 
##     6     2     2     3     2

2. For each value j of this condition, simulate a uniform random variable in the interval \((j,j+1)\). Then make a histogram.

X.5 <- runif(n, min = A.1, max = A.1+1)
hist(X.5, prob = T, breaks =0:20)

3. Compute the sample mean:

mean(X.5)

## [1] 1.499157

3/2

## [1] 1.5

Example 8: Numerical Example

Suppose \((X,Y)\) have the joint pmf

\[ f_{XY}(1,1)= .1, \quad f_{XY}(1,2) = .2, \quad f_{XY}(2,1) = .1, \quad f_{XY}(2,2) = .6, \]

Find the expectations of \(X\) and \(Y\) \[ E(X) = 1.7, \quad E(Y) = 1.8 \]

Find the conditional distributions of \(X\), conditioned on the possible values \(y\) of \(Y\):

\[ \begin{aligned} f_{X|Y}(1|Y=1) &= .5, \quad f_{X|Y}(2|Y=1) = .5\\ f_{X|Y}(1|Y=2) &= .25, \quad f_{X|Y}(2|Y=2) = .75\\ \end{aligned} \]

__Conditional expectations_

\[ E(X|Y=1) = 1.5, E(X|Y=2) = 1.75 \] Formula for \(E(X|Y)\)

\[ E(X|Y) = \begin{cases} 1.5 \quad \text{if} \; Y=1 \\ 1.75 \quad \text{if} \; Y=2 \end{cases} \]

Example 9:

1. Simulation of \((U,V)\)

mydf.0 = data.frame(U = runif(10000), V = runif(10000))
plot(mydf.0,cex = .1)

2. Simulation of \((X,Y)\)

x = sapply(1:10000, function(j){min(mydf.0[j,])})
y = sapply(1:10000, function(j){max(mydf.0[j,])})

# Can also be done with apply()
mydf.1 = data.frame(x =x, y = y)
plot(mydf.1,pch = 1,cex = .1)

3. For the marginal distribution of X, just make a probability histogram

hist(mydf.1$x, prob = T)
abline(a = 2, b = -2, col = 2) 

It shows the straight-line shape that is predicted by theory, $ 2-2x$

4. marginal distribution of Y

hist(mydf.1$y, prob = T)
abline(a = 0, b = 2, col = 2)